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The rate of exchange of C0 2 between the soil and the atmosphere depends on the stability of the organic 
carbon stored in the soil. Recent studies show that carbon stored in the subsoil is characterised by larger 
turnover times than the carbon stored in the topsoil. Consequently, identification of the depth at which high/ 
low amounts of SOC are stored is essential for applying sustainable management of the soil in the light of 
global warming and related threats. This study investigates the depth distribution of SOC in relation to land 
use and soil type based on a large dataset for Flanders (Belgium). Soil type determines the SOC content along 
the entire profile. On the contrary, land use appears to have a strong influence on SOC content in the top 
layers of the profile, but doesn't play a significant role at the bottom of the profile (>1 m depth). SOC content 
near the surface of the profile is remarkably higher in fine (clay) textured soils than in coarse (sand) textured 
soils and tends to increase by increasing soil wetness under sand and silt textured soils. SOC at the bottom of 
the profile increases as well by increasing soil wetness, but only in fine (clay and silt) textured soils. The rate 
of decline of SOC content with depth depends on texture and land use. Under forest this decline is 
remarkably fast, although less so in the more sandy soil types. The overall model predicts the distribution of 
SOC density by depth using land use and soil type information and allows in its integrated form the 
estimation of SOC stocks that can be represented on SOC maps until a reference depth free of choice. 
Applying this model, based on a three-dimensional spatial distribution approach, the total amount of SOC 
stored in Flanders is calculated at 62.20 db 0.72 Mt C for the top 0.3 m and 103.19 ± 1.27 Mt C for the top 1 m. 

© 2009 Elsevier B.V. All rights reserved. 


1. Introduction 

The soil organic carbon (SOC) pool is considered as one of the most 
important reservoirs of the global C-cycle (e.g. Bohn, 1982). This 
reservoir has the potential to act as a major source or sink of 
greenhouse gases due to its large extent and active interaction with 
the atmosphere (Gal et al M 2007). Negative impacts of human activities 
on global warming, soil degradation and desertification can be limited 
by improved soil management, encouraging the enhancement of soil C 
sequestration (Milne et al M 2007). As a consequence, the localisation of 
environments determined by low/high SOC stocks at (sub)national 
level (i.e. the scale where the policy maker acts) is of great importance. 
In addition, organic matter is directly related to soil quality and has a 
major influence on soil structure, water holding and cation exchange 
capacity. Consequently, SOC content is used as an indicator for soil 
quality (e.g. Andrews et al., 2004; Dawson et al M 2007; Pattison et al, 
2008). For example, high soil organic carbon content appears to 
diminish eutrophication or contamination risks of the groundwater or 
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nearby ecosystems due to nitrate (NOJ) or pesticide leaching (e.g. 
Gundersen et al., 1998; Dousset et al., 2004; Poissant et al., 2008). 
Within this context, the vertical distribution of organic matter appears 
to be important, as SOC in the subsoil shows to be strongly related to 
sorption capacity of pesticides and to denitrification capacity of 
leached components (e.g. Lafrance and Banton, 1995; Richards and 
Webster, 1999). These findings underline the importance of mapping 
SOC distribution (with depth) at the regional scale as this information 
can be used as input for models estimating the spatial pattern of 
potential groundwater contamination by pesticides or nitrate leaching 
(e.g. Lindahl et al., 2005; Jarvis et al., 2007). 

The most common approach to map SOC at the regional scale is to 
stratify the study area by land use and/or soil type, to calculate mean 
SOC stocks by land use and/or soil type combination and then attribute 
this value to the corresponding polygons/grid-cells on the map (e.g. 
Kern, 1994; Batjes, 2000; Liebens and VanMolle, 2003; Lettens et al., 
2005). Due to the great spatial heterogeneity of the landscape some of 
these classes are limited in space or not easily accessible and therefore 
risk having few or no SOC data. Consequently, following the traditional 
method no reliable mean SOC value can be predicted for these classes. 
Nevertheless, these classes are most often representing extreme 
environments with extreme (high or low) SOC contents. This implies 
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that these zones have an important contribution to the spatial 
distribution and that the absence of reliable SOC values for these classes 
forms an important source of uncertainty. In more recent studies, a 
model predicting map unit specific SOC values was used in order to 
overcome this problem (e.g. Jones et al M 2004; Cerri et al., 2007; Milne et 
al., 2007; Meersmans et al., 2008). To obtain a more complete and 
detailed spatial distribution Jones et al. (2004) used a pedo-transfer rule, 
defined as a series of ‘if-then’ conditions and with soil, land use and 
climate as input variables, while Meersmans et al. (2008) constructed a 
regression model to predict SOC by land use, texture and drainage 
combination. The latter study showed that no (reliable) predictions 
could be made for large areas of poorly drained valley soils by using the 
traditional mean approach. Consequently, the total area for which 
(reliable) SOC stocks could be calculated increased by 8.1% and the total 
predicted SOC stock increased by 10.1% by applying the model instead of 
the traditional approach. Milne et al. (2007) developed the ‘Global 
Environment Facility Soil Organic Carbon (GEFSOC) modelling system’ 
to map future SOC stock changes at the national and sub-national scale. 
This modelling approach was applied in Jordan (Al-Adamat et al., 2007), 
Brazilian Amazon (Cerri et al., 2007), Kenya (Kamoni et al., 2007) and 
the Indian part of the Indo Gangetic plains (Bhattacharyya et al., 2007). 

Many studies investigated the influence of physical soil properties on 
SOC and/or calculated total SOC storage at regional scale considering a 
fixed soil depth, e.g. top 0.3 or 1.0 m (Feller and Beare, 1997). 
Consequently, these studies do not contain information on the vertical 
distribution of SOC and are missing the depth or 3rd dimension of the 
spatial distribution. Nevertheless, knowledge on the depth distribution 
of soil organic carbon (SOC) is very important to understand how this 
reservoir acts in the global C cycle and in particular to have an insight 
into the soil-atmosphere interaction. In order to investigate the SOC 
distribution with depth, some studies grouped their measurements into 
fixed depth increments (e.g. Jackson et al., 2000; Slobodian et al., 2002; 
Omonode and Vyn, 2006), while others fitted continuous functions 
through the data. Thereby, the exponential decline is most frequently 
used (e.g. Hilinski 2001 ; Sleutel et al., 2003 ). The factors determining the 
SOC distribution with depth are various. Different studies showed that 
root biomass declines faster with depth than SOC (e.g. Gill et al., 1999; 
Jobbagy and Jackson, 2000). This should be the consequence of a 
downwards C transport within soil profile or lower decomposition rates 
of organic material in deeper layers. On the one hand, seepage water 
transport, texture related soil permeability, earthworm bioturbation, 
leaching and vertical soil mixing by organisms were identified in 
previous studies as the most important factors determining the vertical 
transport of SOC (Jobbagy and Jackson, 2000; Don et al., 2007). On the 
other hand, radiocarbon studies illustrate an increase in turnover time 
and stability of SOC with increasing depth (e.g.Trumbore, 2000; Baisden 
and Parfit, 2007). Fontaine et al. (2007) considered the lack of energy 
supply of microbes in deeper soil layers, under the form of fresh organic 
matter as the main explanatory factor for this phenomenon. Thus, 
defining the depth at which important SOC densities are present 
indicates the stability of the carbon stored in the soil. In other words, SOC 
information by depth provides information on the likelihood of the soil 
to exchange C0 2 with the atmosphere. 

This study aims to construct an empirical model describing depth 
distribution of soil organic carbon (SOC) in relation to texture, 
drainage and land use for Flanders (Belgium) in order to analyse the 
spatial distribution of SOC density in three dimensions and to produce 
SOC stock maps until a reference depth to be chosen by the user. 

2. Material and methods 

2.2. Study area 

Flanders is situated in the northern part of Belgium, a country in 
the north-west of Europe, and covers an area of 15,521 km 2 . Given the 
rather small size and the absence of important variation in altitude, 


this region is characterised by uniform climate conditions, with an 
annual precipitation amount of 700 to 800 mm and a mean 
temperature of 9 to 10 °C. Nevertheless, one can find a variation in 
soil types in Flanders: fertile loess soils in the south (Luvisol), sand 
textured soils in the north (Podzols) and wet clay rich soils (Fluvisols) 
in the alluvial and coastal plains in the west. As the region is densely 
populated, the area is dominated by agricultural and urban land use 
types. 

2.2. SOC data 


The model has been developed based on the dataset ‘Aardewerk’, 
containing almost 7000 profiles pits sampled throughout Flanders 
during the Belgian National Soil Survey (1947-1974). 75% of the data 
was used to calibrate the model. The other 25% of the data were used for 
validation. The profiles in the dataset were provided with information 
on land use, texture and drainage class, classified according to the 
Belgian soil classification system (Fig. 1 and Table 1 ). No information on 
land use history of the sites was reported. In most of these profiles, the 
soil is described until a depth of 1-1.5 m, and contains 5 to 6 horizons. 
For each horizon the depth, the percentages of sand, silt, clay, as well as 
the SOC concentration (g C kg -1 ), determined by the classic bichromate 
method of Walkley and Black (1934) are given (Van Orshoven et al., 
1988). In order to compensate the incomplete oxidation of carbon 
characteristic for the Walkley and Black method, a correction factor of 
1.33 is applied. SOC mass densities (kg m~ 3 ) are obtained by 
multiplying the SOC mass concentration (g C kg - ^ by soil bulk density 
(kgrrT 3 ) (Eq. (1)): 


SOC = p s * 


C 

Tooo 


where: 

SOC SOC mass density of the sample (kg C m~ 3 ) 
ps bulk density of the soil (kg m -3 ) 

C SOC concentration of the sample (g C kg -1 ). 



As bulk density information is not provided by ‘Aardewerk’, a 
pedo-transfer function (PTF) is used to predict this soil property for 
each sample. The use of a PTF risks to over/underestimates real bulk 
density and creates a source of uncertainty in calculating SOC mass 


U Heavy clay 
EClay 
A Silt loam 

L (heavy) Sandy loam 
P Light sandy loam 
S Loamy sand 
Z Sand 


100 



Silt content (%) 
(2pm - 50pm) 


Fig. 1. Belgian soil texture classification triangle (after Ameryckx et al., 1995). 
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Table 1 

Belgian drainage classification system as a function of the soil texture (heavy clay (U), 
clay (E), sandy loam (L), silt loam (A), light sandy loam (P), loamy sand (S), sand (Z); 
see Fig. 1 for the Belgian texture triangle) after Ameryckx et al. (1995). 


Drainage 

class 

Definition 


Depth oxidation 
horizon (m) 

(min. depth 
water table) 

Depth reduction 
horizon (m) 
(max. depth 
water table) 


Belgian system 3 

USDA 

drainage class 5 

Text. A, 

L, E, U 

Text. 

4 S, P 


a 

Very dry 

Excessively 

- 

- 

- 

b 

Dry 

Well 

>1.2 

0.9-1.2 

- 

c 

Moderately dry 

Moderately- 

well 

0.8-1.2 

0.6-0.9 

— 

d 

Moderately 

Wet 

Moderately 

0.5-0.8 

0.4-0.6 

— 

h 

Wet 

Moderately- 

poorly 

0.3-0.5 

0.2-0.4 

— 

i 

Wet 

Poorly 

0-0.3 

0-0.2 

- 

e 

Wet with 
reduction 
horizon 

Moderately- 

poorly 

0.3-0.5 

0.2-0.4 

>0.8 

f 

Very wet 
with reduction 
horizon 

Poorly 

0-0.3 

0-0.2 

0.4-0.8 

g 

Extremely wet 

Very poorly 

0 

0 

<0.4 


a Ameryckx JB, Verheye W, Vermeire R (1995) Bodemkunde, Gent. 

5 Soil Survey Staff (1951) Soil Survey Manual, United States Department of 
Agriculture, U.S. Dept. Agricultural Handbook, 18, Washington, D.C. 


SOC database correspond to the one of digital soil map, which is 
important to clearly identify the influence of different soil properties 
on SOC distribution with depth. 

2.4. SOC depth distribution 

SOC depth distribution in the tillage layer differs from the deeper 
layers. The SOC density remains constant until the tillage depth (td) is 
reached. At greater depth, SOC density shows an exponential decline 
with depth (Eq. (3)): 

z < td : SOC(z) = SOC surf 

z > td : SOC(z) = A • e“' (z_td) + SOC 1 ’ 

where: 

z depth (m) 

td tillage depth (m) 

SOC(z) SOC mass density at depth z (kg C m~ 3 ) 

SOC surf SOC mass density (kg C m -3 ) at the surface 

SOCoo SOC mass density (kg C m~ 3 ) at the bottom of the soil profile. 

a is a constant which determines the shape of the exponential part 
of the curve. 

A — SOC surf — SOC. (4) 


density. Nevertheless, Boucneau et al. (1998) compared different 
pedo-transfer functions (PTF's) predicting soil bulk density for 
Flemish soils. They showed that, of all PTF's (applicable in this 
study) the general PTF of Manrique and Jones (1991), has the best 
correlation between observed and predicted bulk density values. 
Consequently this function, with a minimal value of 0.3 g cm -3 (i.e. 
value corresponding to peat soils (Lettens et al., 2004)), is used in this 
study to predict bulk density (g cm -3 ) for each sample (Eq. (2)). 

Ps = 1.66-0.318*^ (2) 


2.3. Land use and soil type maps 

The study area was stratified by an overlay of the digital land use map 
of Flanders with the soil map of Belgium. The digital land use map of 
Flanders, with a resolution of 15 m, was derived from Landsat7-ETM+ 
images of the year 2001. Classification of these images combined with 
external road and waterway information resulted into 18 land use 
classes (OC GIS Vlaanderen, 2002). For the present study, four 
aggregated classes were extracted from the digital land use map of 
Flanders: i.e. forest, grassland, cropland and heath. 

The digital soil map of Belgium is derived from the National Soil 
Survey (1947-1974) (OC GIS Vlaanderen, 2001). From the soil codes 
on this map information concerning drainage and texture class could 
be extracted. Fig. 1 illustrates that Belgian soil texture classes were 
defined according to clay, silt and sand content. Whereas Table 1 
shows that the Belgian soil drainage class system is based on the depth 
of occurrence of oxidation and reduction properties in the soil profile. 
These depths correspond to the minimum and maximum depth of the 
ground water during the year, or the position of the ground water 
table in winter and summer, respectively (Ameryckx et al., 1995). 
Consequently, Belgian texture and drainage classification systems 
allow production of numeric maps presenting variables related to 
these properties, i.e. regarding texture: clay, silt, and sand content as 
well as geometric mean particle size (Dg) and considering drainage: 
maximal and minimal depth of the ground water table. As drainage 
information was not updated for current situation, inventories of the 


When the model is applied to permanent grassland, forest or 
heath, Eq. (3) can be simplified to an exponential relation: 

SOC(z) = A • exp“' z + SOC„. (5) 

In order to obtain a realistic tillage depth (td), different model 
simulations were performed for cropland. For each simulation the 
data was divided in two groups: i.e. one group with samples near the 
surface to fit the constant part of the model and another group, 
containing all deeper samples, to fit the exponential decline at greater 
depth. In each simulation the amount of samples in the first group 
increased and in the second group declined one at the time. The 
simulation with the best fit was retained. From the selected 
simulations the tillage depth (td) is determined by the intersection 
of the constant and exponential parts of the curve (Eq. (6)). More 
detailed explanation and a graphical example of this method can be 
found in Meersmans et al. (in press). 

\n ( ( S0Csurf ~ S0C ~A 

^ — V A ' ) fC 2\ 


where 

A! extrapolated SOC(O) of the exponential part of the model — 

SOCco. 

2.5. SOC mass 

To calculate SOC mass per surface unit (kg C m~ 2 ) until a certain 
reference depth (rd) Eq. (3) is integrated and results in Eq. (7): 

A ■ (e“'(rd-td) _ l) 

SOC mass = SOC surf • td + —---- L + SOC. • (rd - td). (7) 

Fig. 2 presents the methodological flowchart of this study. After 
grouping the data by land use-soil type class, the model (Eq. (3)) was 
fitted through the SOC depth data for all land use-soil type classes 
with more than 30 samples (Fig. 3). This resulted into a set of SOCoo, A, 
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SOC database Aardewerk 



ANOVA Depth distribution (eq. 3) 

(by depth increment) (by land use - soil type subdataset) 


\ 

\ 

ii 

Empirical Models (PTF's) 

parameters depth distribution variables describing 

part lb (A, a, SOCoo, td) ^ environment 

(land use, soil type) 


Substitution of empirical models (PTF's part II) in eq. (3) 

III 

One physical 3D SOC model 

SOC - f(depth, land use, soil type) 




Substitution of empirical models (PTF's part II) in eq. (7) 

iv_\j/_ 

SOC map + total stock (Flanders) 

Reference depth free of choose 
(example 0,3 and 1 meter) 


td and a parameter values for each environment, which allows 
construction of a PTF for each parameter of the depth distribution 
function (Eq. (3), A , SOCo, a, td) by expressing the parameter as a 
function of land use, drainage and soil texture variables. To obtain one 
physical model describing depth distributions of SOC as a function of 
land use and soil type, the resulting pedo-transfer functions, are 
substituted in Eq. (3). Furthermore, substitution of these four PTF’s in 
Eq. (7), i.e. the integrated form of the depth distribution model, 
combined with digital land use and soil type maps of the study area 
allows calculation of the total amount of carbon stored in the study 
area. By integrating the depth distribution function until certain depth 
of interest, the reference depth will be chosen by the user. 


2.6. Software 

For the modelling part of this research Matlab version 7.3 
(MathWorks, Natick, Massachusetts, USA) was used. Statistical 
analyses were carried out in SPSS 11.0 (SPSS Inc., Chicago, Illinois, 
USA). 


3. Results and discussion 

3.1. ANOVA by depth increments 

As expected, the results of ANOVA tests after grouping the data by 
20 cm depth increments indicate that the influence of land use on SOC 
density decreases with depth, as land use has significant effect 
(p = 0.001) on SOC in the top layers, but not in the subsoil (>80 cm; 
Table 2). Texture and drainage, on the other hand, have a significant 
(p = 0.001) effect on the SOC content throughout the profile. More¬ 
over, the interaction between texture and drainage is significant 
(p = 0.001) for all depths, which indicates that the effect of texture on 


Fig. 2. Flowchart of the methodology. 


Texture class 


Land Use 


Cropland 


Grassland 


Forest 



0 10 20 30 40 50 60 70 0 10 20 30 40 50 60 70 0 10 20 30 40 50 60 70 80 


SOC mass density (kg m 3 ) 


Fig. 3. Depth distribution model (Eq. (3)) applied on different land use-drainage class combinations within clay (E), silt-loam (A) and sand (Z) textured soils (n>30). A complete 
figure whereby the model is applied on all land use-soil type combinations of database Aardewerk (n>30) can be obtained from following website (i.e. E-Figure 1): www. 
jeroenmeersmans.be/research/publications/depth3Dmod. 
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Table 2 

ANOVA results, testing the influence of land use (LU), texture (TXT), drainage (DR) and 
interactions on the SOC stocks by depth increments of 0.2 m at 3 different levels of 
significance: 0.01<p<0.05 (*), 0.01 <p<0.001 (**) andp<0.001 (***) (N indicates the 
number of samples). 



SOC density depends on the drainage status (and visa versa) for all 
depths. These trends are in agreement with the ones of Jobbagy and 
Jackson (2000) who showed that the SOC content at the surface (0- 
20 cm) appears to be dominated by land use and climate, i.e. strong 
positive relation with precipitation and negative relation with 
temperature. While in contrast, SOC shows the highest correlation 
with clay content in deeper layers (1-3 m). 

3.2. Depth distribution model 

The depth distribution model (Eq. (3)) was applied for all land 
use-soil type combinations with at least 30 samples (E-Figure 1; 
Meersmans, 2008). Each graph shows the SOC depth distributions of 
different drainage class groups within one land use-texture class 
combination. A selection of some land use (i.e. cropland, grassland, 
forest)-texture classes (i.e. clay (E), silt loam (A) and sand (Z)) of E- 
Figure 1 is given in Fig. 3. From this figure one can observe 
immediately the general known trends to have higher SOC mass 
density values under more humid and/or fine textured soils and that 
under cropland tillage diminish considerably SOC content in the top 
layer (e.g. Tan et al., 2004; Zinn et al., 2005). By investigating the 
depth distribution instead of SOC content until fixed depth, we were 
able to identify the depth of specific relationships between these soil 
variables and the SOC content. Therefore these environmental 
variables will be discussed in detail for each parameter of the general 
SOC distribution model (Eq. (3)). 

3.2.1. PTF A-parameter 

The results show that topsoil SOC is highest under forest and lowest 
under cropland in almost all land use-soil type combinations (Fig. 3, E- 
Figure 1). This is in accordance with the literature and explained by the 
fact that forests have a relative higher above ground allocation of organic 
material as compared to other land uses. This results in a higher input of 
SOC at the surface due to litter fall under forest (Wang et al., 2004). 
Moreover, recent research pointed that difference in land management 
practises could affect seriously the distribution of SOC with depth. Under 
cropland, tillage declines the SOC stock in the top layer by exposing 
periodically the organic matter that is physically protected in micro¬ 
aggregates to biodegradation (Balesdent et al., 2000). Under grassland, 
the SOC mass densities near the surface tend to increase with decreasing 
particle size (i.e. higher clay or silt content). This is in agreement with 
other studies reporting a positive linear relation between clay or clay + 
silt content and SOC content until fixed reference depth (e.g. Zinn et al., 
2005). A positive relation between soil wetness and SOC mass density 
near the surface can be observed in the sand and silt textured grassland 
soils (A, L, P, S, and Z). This trend was not found under clay textured soils 
(E, U) (Fig. 3, E-Figure 1). Davidson (1995) underlines as well the 
importance of drainage status on SOC by showing that SOC content 
increases with soil wetness. The relation between SOC mass near the 
surface and soil type is different for cropland and grassland. Under 
cropland the heavier soils (E, U, L) show a positive relation between SOC 
mass density near the surface and soil wetness, whereas this trend 


cannot be observed in sand and silt loam soils (A, Z) (Fig. 3, E-Figure 1). 
As the A-parameter is strongly correlated with the SOC mass density 
near the surface (i.e. correlation coefficient between A and SOC surf of 
0.851 for cropland, 0.966 for grassland, 0.999 for forest and 0.999 for 
heath, data not shown) the extracted topsoil SOC trend information is 
very useful to construct the PTF describing the A-parameter (Eq. (3)) as 
function of land use and soil type. The PTF contains 3 different land use 
specific equations (Eq. (8)), because of the strong land use dependency 
of the SOC mass density near the surface and their relation with soil 
type: 

^grassland = H 2 O s * (a*silt 2 + f>*silt) + c*H 2 O s *sand 
+ d*H 2 O w + e*clay +/ 

^cropland = g*H 2 O s *(clay + silt) + ft*sand 2 + i*sand + j*clay + k 

^forest/health = **H 2 O s + ™*Dg + H 

( 8 ) 

where: 

H 2 O s height of the deepest ground water table position (summer) 
above 1.5 m depth (m) 

H 2 O w height of the shallowest ground water table position 
(winter) above 1.5 m depth (m) 
sand sand content (%) 

clay clay content (%) 

silt silt content (%) 

Dg geometric mean particle size diameter (mm). 

The corresponding R 2 value of the model is 0.81 for all environ¬ 
ments, with land use specific values of 0.60 for grassland, 0.61 for 
cropland and 0.50 for forest/heath. Fig. 4 compares for every land 
use-soil type combination the initially calculated A-parameter value 
by applying Eq. (3) with the predicted values based on the PTF 
(Eq. (8)). This figure illustrates that the A-parameter values under 
cropland are remarkably lower than under grassland, forest/heath. 

3.2.2. PTF SOCoo-parameter 

One can state that for sand and sand-loam (L, P, S, Z) soils for each 
drainage-land use combination the SOC mass density converges to 
zero near the bottom of the profile (SOC,), while for clay and silt loam 
soils, this value is remarkably higher and tends to increase with 
increasing soil wetness (Fig. 3, E-Figure 1). Within the same soil type 
no important differences in SOCoo could be observed between the 



Fig. 4. Calculated (Eq. (3)) versus PTF predicted (Eq. (8)) A parameter (kg m 3 ) for all 
land use-texture-drainage class combinations (Fig. 3). 
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Fig. 7. Mean and standard deviation of calculated (Eq. (3)) versus PTF predicted (Eq. (11)) 
tillage depth by texture class (i.e. U: heavy clay, E: clay, A: silt loam, L: (heavy) sandy loam, 
P: light sandy loam, S: loamy sand, Z: sand). 


T I 

o 

n 



• calculated 
a predicted 

ii i i i i 

U E A L P S Z 

texture class 


Fig. 5. Calculated (Eq. (3)) versus PTF predicted (Eq. (9)) SOCoo parameter (kg m 3 ) for 
all land use-texture-drainage class combinations (Fig. 3). 

different land uses. This corresponds to ANOVA analysis discussed 
earlier (Table 2). Moreover, the results lies in the line of the finding's 
reported by Jobbagy and Jackson (2000) showing that clay content 
controls SOC content in deeper layers. Consequently, a single 
relationship between SOCoo and soil type (i.e. texture and drainage) 
could be found. This resulted in one model, predicting SOCoo, for all 
land uses (Eq. (9)). Moreover, the important interaction between 
texture (clay and silt content) and drainage (H 2 O s ) on SOC stock at 
these depths is highlighted in the equation. 

SOC^ = a*H 2 O s *clay ± b*H 2 O s *silt ± c*clay ± d*silt (9) 

The model (Eq. (9)) nicely describes the observed pattern with an 
R 2 value of 0.76 (Fig. 5). 

3.2.3. PTF a-parameter 

The rate of decline of SOC with depth in the exponential part of the 
depth distribution model (Eq. (3)) is represented by the ct-parameter. 
As expected, the results show that forested soils are characterised by 
high SOC stocks near the surface and a strong decline of SOC with depth, 
while grassland soils have lower SOC stocks in top layers but higher 
stocks in deeper layers as a result of a weaker decline of SOC with depth 


sand % 

100 90 80 70 60 50 40 30 20 10 0 



Fig. 6. Mean and standard deviation of calculated (Eq. (3)) and PTF predicted (line) 
(Eq. (10)) a parameter by sand content (%). 


(Fig. 3, E-Figure 1). This can be explained by the difference in below 
versus above ground allocation of living organic matter between both 
land uses. Mokany et al. (2006) measured for temperate climate root to 
shoot ratios of 0.20 ±0.03-0.46 ±0.06 under forest and 4.22 ±0.52 
under grassland. As under forest, almost all living organic material is 
stored above ground in the trees, the input of SOC mainly takes place at 
the surface due to fall of woody and relative slow decomposable debris. 
By way of contrast, under grassland the majority of the living organic 
material is stored below ground in the well developed root systems, 
which cause an important input of SOC in the subsoil from root turnover 
(Jobbagy and Jackson, 2000; Wang et al., 2004). Under sand textured 
soils the decline of SOC with depth in the subsoil is less pronounced than 
in heavier soils (Fig. 3, E-figure 1). This is counter intuitive as coarse 
textured soils are characterised by larger pores than fine textured soils 
and so have a higher oxidation capacity of SOC, resulting in a shorter 


Table 3 

Parameter values and 95% confidence interval of PTF's describing A, a, SOCoo and td as a 
function of land use and soil type (Eqs. (8)—(11)). 


Parameter depth 
distribution model 
(Eq. (3)) 

Land use 

Parameter Value 

PTF 

95% confidence interval 

A (kg m~ 3 ) 

Grassland 

a 

0.0056 

0.0020 

0.0091 



b 

-0.3269 

-0.5945 

-0.0592 



c 

0.1454 

0.0766 

0.2143 



d 

11.7340 

9.1525 

14.3154 



e 

0.4135 

0.3431 

0.4839 



f 

21.7186 

19.5789 

23.8584 


Cropland 

g 

0.1522 

0.0823 

0.2221 



h 

0.0013 

0.0010 

0.0016 



i 

-0.0665 

-0.0925 

-0.0406 



j 

0.1304 

0.0916 

0.1691 



k 

14.315 

13.8036 

14.8263 


Forest 

l 

42.3713 

34.7906 

49.9511 



m 

-29.1021 

-33.5331 

-24.6710 



n 

53.7206 

51.4357 

56.0054 

a 

Grassland 

a 

0.0155 

0.0173 

0.0138 



b 

-4.8747 

-4.7674 

-4.9819 


Cropland 

c 

-0.0003 

-0.0003 

-0.0004 



d 

0.0448 

0.0504 

0.0393 



e 

-4.8863 

-4.7907 

-4.9819 


Forest 

f 

0.0337 

0.0367 

0.0308 



g 

- 6.9729 

-6.7565 

- 7.1893 

SOCo (kg m~ 3 ) 

All 

a 

0.0794 

0.0112 

0.1476 



b 

0.1226 

0.1033 

0.1419 



c 

0.1216 

0.1121 

0.1311 



d 

0.0162 

0.0147 

0.0177 

td (m) 

Cropland 

a 

-0.0010 

-0.0017 

-0.0003 



b 

0.1018 

0.0709 

0.1327 



c 

0.1219 

0.1103 

0.1335 
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Relative error SOC mass density (%) 


Fig. 8. Frequency distribution of the relative error on predicted SOC mass densities by applying the overall depth distribution model to different soil depths, land uses and soil type 
conditions. 


turnover time. This implies that the slower decline of SOC with depth in 
sandy soils can not be the consequence of more stable conditions in 
deeper layers. An alternative explanation would be that vertical 
transport of SOC must be an important process for C accumulation at 
these depths. Don et al. (2007) stated that this process is stronger in 
profiles characterised by higher soil permeability. The results of our 
study confirm this, as the SOC decline with depth tends to slow down 
with increasing sand content, which suggests a stronger vertical SOC 
transport in soils characterised by larger pores or a higher permeability. 
The ai-parameter is expressed in this study as a function of sand content 
in three different land use specific equations (Eq. (10)). Under grassland 
and forest the relation between the ai-parameter and sand content could 
be described by linear functions. Under cropland the expression is more 
complex, as a quadratic function better describes the system. 

^-grassland fl*Sand + b 

^cropland = c*sand 2 + d*sand + e (10) 

^forest = /* sand + g 

Fig. 6 illustrates the ct-parameter by sand content and by land use 
class. The three land use specific models appear to converge to the 
same value (— 3.5 to — 4) for a sand content of 100%. As near all heath 
soils appear on sandy soils in Flanders (E-Figure 1), this land use could 
only be studied under Z textured sand (90% of sand). They are, with a 
mean ai-parameter value of — 5.7 ± 0.3, characterised by a very strong 
decline of SOC with depth. Despite the clear trends one can observe in 
Fig. 6, the proposed PTF has a rather low R 2 of 0.42, with R 2 values for 
the land use specific models of 0.23 for cropland, 0.20 for grassland, 
and 0.53 for forest. 

3.2.4 PTF tillage depth 

The calculated tillage depth (td) depends on texture. The mean 
tillage depth by texture class varies between 8 cm for heavy clay (U) 
and 19 cm for sand soils (Z; Fig. 7). This trend was expected and can be 
easily explained by the fact that the traction required for ploughing in 
fine textured soils (E, U, A) is much larger than in coarse textured soils 
(S, Z) and animals or tractors used in the 1960s were probably not 
powerful enough to plough the same depth in fine as in coarse 
textured soils. Increased mechanisation during the last 50 years allows 
deeper ploughing of the soil. As a consequence much deeper tillage 
should be expected nowadays, foremost in the fine textured soils. Van 
Meirvenne et al. (1996) recorded for northwest Belgium an increase in 
average tillage depth of almost 0.10 m between 1960 and 1990. 
Meersmans et al. (in press) showed that this increase in tillage depth 


in Belgian cropland soils is more pronounced under fine textured than 
coarse textured soils. They calculated an increase of 0.079 ± 0.031 m in 
silt loam (A) soils and 0.024 ± 0.050 m in sand (Z) soils for the period 
1960-2006. Tillage depth is strongly correlated with Dg and clay 
content. This is shown by a correlation coefficient (r) of —0.65 and 
0.69 respectively (data not shown). As the correlation between clay 
content and Dg was not too strong (<0.6) a linear combination of both 
texture variables is used to predict tillage depth (Eq. (11)). 

td = a*clay + b* Dg + c (11) 

The corresponding R 2 value of the model is 0.57. Nevertheless, Fig. 7 
illustrates that the predicted tillage depths using the model (Eq. (11)) 
are close to the observed mean tillage depths for each texture class. 

3.3 Overall (land use-soil type) SOC depth distribution model 

The substitution of the A, SOCo, td and a parameter PTFs (Eqs. (8)- 
(11 )) in the depth distribution model (Eq. (3)) and calculation of SOC surf 



Fig. 9. Validation of the overall depth distribution model by comparing mean observed 
and predicted SOC mass densities (kg m~ 3 ) of the validation sub-dataset. Each point 
represents a land use-soil type-depth increment (0.1 m interval) combination with at 
least 5 samples. 
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using Eq. (4), allows calculation of SOC density as a function of depth, 
land use and soil type. In other words the overall model allows 
generation of SOC distributions with depth for a given land use-soil type 
combination. The parameter values as well as their 95% confidence 
limits of this overall depth distribution model are given in Table 3. All 
parameters are significant at p< 0.05 level. 

3.3.1. Model uncertainty 

In order to calculate the error on the estimated SOC mass density 
(kg m -3 ) due to parameter uncertainty (Table 3) of this overall depth 
distribution model the Monte Carlo simulation technique was applied. 
Hereby, the model predicts the SOC mass density (kg m~ 3 ) for each 
simulation by a set of parameter values randomly chosen from normal 
probability distributions with mean and standard deviation equal to the 
estimated value and 66.7% uncertainty (stdev) of the corresponding 
parameters. From all predictions, the standard deviation is used as error 
for the estimated SOC mass density (kg m~ 3 ). The error on the 
estimated SOC mass density (kg m~ 3 ) strongly depends on land use, 
texture, drainage and soil depth. Regarding the top 1 m, the relative 
value of this error ranges between 1.86% (for dry silt loam cropland soils 
in the tillage layer) and 28.65% (for extremely wet grassland at the soil 
surface). The majority (>80%) of the relative errors are between 2% and 
7% (Fig. 8). 

3.3.2. Model validation 

To evaluate the quality of this overall (land use-soil type) depth 
distribution model we compared predicted and mean observed SOC 
mass density (kg m -3 ) by land use-soil type-depth increment (0.1 m 
interval) combination with at least 5 samples from the validation set 
(Fig. 9). The validation indicates that the model is able to predict 
reliable values for most environments, with an associated RMSE of 
3.27 kg SOC m~ 3 and a ratio of performance to deviation (RPD), i.e. 


the ratio of the standard deviation to RMSE, of 2.89. Nevertheless, one 
can observe an underprediction of the model for some land use-soil 
type-depth increment combinations (Fig. 9). Two different environ¬ 
ments for these outliers can be distinguished, i.e. on the one hand, 
subsoil horizons most likely containing peat and situated in 
moderately to very poorly drained soils, and on the other hand, Bh 
horizons enriched with organic matter (depth around 0.3-0.4 m) in 
well to moderately drained sandy podzols. 

3.4. Model application — case study: SOC inventory Belgium 

The estimated spatial distributions of SOC stocks in Flanders until 
reference depths of 0.3 and 1 m are calculated by combining the 
integrated form of the overall depth distribution model (Eq. (7)) with 
the PTF expressions (Eqs. (8)—(11 )) and digital land use and soil type 
maps of the region (Figs. 10 and 11). Besides the importance of these 
maps regarding the understanding of current climate change process, 
they can be useful as input in estimating other important environmental 
threats, like potential groundwater contamination by pesticides or 
nitrate. As expected, the pattern of valleys and river systems can be 
recognised on both maps. The reduced oxidation capacity of SOC under 
anaerobic conditions results in remarkably higher carbon contents in the 
poorly drained valley and depression soils compared to the dry soils on 
the plateau (Meersmans et al., 2008). Moreover the clay rich soils in the 
coastal planes are characterised by high SOC densities. In these soils, 
organic matter is physically protected against degradation by microbial 
attack through sorption of organic matter to clay minerals and enclosure 
within the soil aggregates (Razafimbelo et al., 2008). While the map 
representing the SOC stock until 1 m depth (Fig. 11) is entirely 
dominated by soil type, the map showing SOC until a reference depth 
of 0.3 m (Fig. 10) also reflects the land use pattern. On this map, forested 
areas in the north-eastern part of the region as well as south to south- 
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Fig. 10. SOC distribution map for Flanders (Belgium) for the top 0.3 m (kg m 2 ). 
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Fig. 11. SOC distribution map for Flanders (Belgium) for the top 1.0 m (kg m 2 ). 


east of Brussels can be easily recognised through their higher SOC 
densities in the topsoil. The white spots on the maps correspond to 
major cities (e.g. Antwerp and Ghent) or military areas where no soil 
type data is available. Consequently, only 77.1% of the total area of 
Flanders (i.e. 10,428 km 2 ) could be taken into account to calculate the 
total SOC stock for Flanders. The total SOC stocks for the corresponding 
area are estimated at 103.19 d= 1.27 Mt C until a reference depth of 1 m 
and 62.20 =b 0.72 Mt C until a reference depth of 0.3 m. 

These results are in agreement with the ones obtained by Meers¬ 
mans et al. (2008) using a multiple regression model for the top 1 m soil. 
These authors also pointed to the relation between the fluvial pattern 
and the spatial distribution of SOC in Flanders. Moreover the mean SOC 
mass densities by land use predicted following the present 3D model 
approach are comparable to the mean values for 1960 from studies 
following the common approach (i.e. calculating mean SOC stocks by 
landscape unit). The SOC mass density under grassland in this study 
(11.5 kg C m -2 ) fits well in the range of mean values reported by other 
studies, e.g. Liebens and VanMolle (2003) (12.9 kg C m~ 2 ) and Lettens 
et al. (2004) (11.7 kg m -2 ), while the SOC mass density under cropland 
in this study (7.9 kg C m -2 ) is situated in between the values obtained 
by Liebens and VanMolle (2003) (7.4 kg Cm -2 ) and Lettens etal. (2004) 
(8.9 kg C m~ 2 ). All values are close to the overall mean SOC mass density 
value under cropland (7.8 kg C m -2 ) reported for Flanders in 1990 by 
Sleutel et al. (2003). 

4. Conclusions 

The results of this study demonstrate that the SOC stock near the 
surface is determined by land use and soil type, while SOC near the 
bottom of the soil profile only depends on soil type (i.e. texture and 
drainage). Moreover, the interaction between texture and drainage has 
an important influence on the distribution of SOC with depth. In general, 


under grassland, increasing soil wetness results in higher SOC densities 
near the surface under sand and silt textured soils (Z, S, P, L, A) and into 
higher SOC densities at the bottom of the profile under silt and clay 
textured soils (U, E, A). Furthermore, under cropland higher clay content 
results in remarkable higher SOC contents throughout the profile in the 
poorly drained soils. The shape of the curve depends on the land use 
and texture. Under cropland, SOC remains constant in the tillage layer 
and is remarkably lower than at the same depths under grassland or 
forest. Tillage depth appears to be shallower in fine textured clay soils 
and tends to increase towards the sandy textured soil types, because in 
the latter ploughing was less limited by the traction in the 1960s 
(period of the database). Lower root to shoot ratios under forest 
compared to other land uses, resulting in a relative higher input of 
organic carbon by leaf fall (than at subsurface by roots), explain the 
more pronounced SOC declines with depth under forest than under 
grassland or cropland. Moreover, the SOC decline with depth is reduced 
when sand contents increase, suggesting that vertical transport of SOC 
is promoted in the more permeable sandy soils. The resulting SOC map 
for Flanders shows a strong correlation with the pattern of river valleys 
in the study area, as the poorly drained soils are characterised by high 
SOC contents. Furthermore, the clay rich soils in the coastal planes are 
also characterised by high SOC densities. The total SOC stock in the area 
is estimated at 62.20 ± 0.72 Mt C and 103.19 ± 1.27 Mt C for top 0.3 and 
1 m of the soil, respectively. A future perspective can be to extend the 
present model to account for land use history or changed drainage 
conditions. 
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Appendix A. Supplementary data 

Supplementary data associated with this article can be found, in 
the online version, at doi:10.1016/j.geoderma.2009.05.015. 
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